function van

global lam

figure(1)
clf
% get(gcf)
set(gcf,'Position', [25 1028 560 309])
co = [0 0 1;
    0 0.5 0;
    1 0 0;
    0 0.75 0.75;
    0.75 0 0.75;
    0.75 0.75 0;
    0.25 0.25 0.25];
set(groot,'defaultAxesColorOrder',co)
subaxis(1,1,1,1,'MT',0.01,'MB',0.2,'MR',0.01,'ML',0.05,'P',0.04)

lam=10;

hold on
%axis equal
n=100;
th=linspace(0,2*pi,n);
for i=1:n
    ys(i)=2*cos(th(i));
    vs(i)=2*sin(th(i));
end

tmax=20*pi
y10=3; y20=3;
y0=[y10 y20];
%  calculate solution using a MATLAB routine
options = odeset('RelTol', 1e-12,'AbsTol',1e-12);
[t,y] = ode45(@rhs,[0 tmax],y0,options);

tmax=20*pi
y10=0.5; y20=0;
y0=[y10 y20];
%  calculate solution using a MATLAB routine
[t2,y2] = ode45(@rhs,[0 tmax],y0,options);

% determine limit cycle
tmax=100*pi
y10=1; y20=0;
y0=[y10 y20];
%  calculate solution using a MATLAB routine
[tl,yl] = ode45(@rhs,[0 tmax],y0,options);
k=[];
for i=2:length(tl)
    if yl(i,2)<0 & yl(i-1,2)>0
        k=[k i];
    end
end
nk=length(k);

% find period
for i=2:nk
    Tl=tl(k(i))-tl(k(i-1));
end


hold on
plot(y(:,1),y(:,2),'r','LineWidth',1)
plot(y(1,1),y(1,2),'.r','MarkerSize',20)

plot(y2(:,1),y2(:,2),'b','LineWidth',1)
plot(y2(1,1),y2(1,2),'.b','MarkerSize',20)

%plot(ys,vs,'-.k','LineWidth',1.5)
plot(yl(k(nk-1):k(nk),1),yl(k(nk-1):k(nk),2),'--k','LineWidth',1.5)

%text(-3,3.2,'$\lambda =  0.1$','Interpreter','LaTex','FontSize',16)
i=200;  ii=i+2;
%arrowhead([y(i,1) y(ii,1)],[y(i,2) y(ii,2)],'r',[2.8 1.4],2);

i=900;  ii=i+150;
%arrowhead([y2(i,1) y2(ii,1)],[y2(i,2) y2(ii,2)],'b',[2.1 1.5],2);

%axis([-3.5 4 -4 4])

% commands to label each axes
xlabel('y-axis','FontSize',12,'FontWeight','bold')
ylabel('v-axis','FontSize',12,'FontWeight','bold')
grid on

% have MATLAB use certain plot options (all are optional)
box on
% Set the fontsize to 14 for the plot
set(gca,'FontSize',14)

hold off


%exportgraphics(gcf,'/Users/mark/Desktop/fig616.pdf')

figure(2)
clf
% get(gcf)
set(gcf,'Position', [26 524 560 420])
co = [0 0 1;
    0 0.5 0;
    1 0 0;
    0 0.75 0.75;
    0.75 0 0.75;
    0.75 0.75 0;
    0.25 0.25 0.25];
set(groot,'defaultAxesColorOrder',co)
subaxis(1,1,1,1,'MT',0.01,'MB',0.2,'MR',0.01,'ML',0.05,'P',0.04)

hold on
plot(t2,y2(:,1),'b')
plot(t2,y2(:,2),'r')



%  define f1(t,y)
function dy=rhs(t,y)
global lam
dy=zeros(2,1);
dy(1) = y(2);
dy(2) = lam*(1-y(1)^2)*y(2)-y(1);









